#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <sys/time.h>
#include <gmp.h>

#define BIT_PRECISION 33219281

int main(void)
{
	struct timeval start, end; // Variaveis para guardar o tempo de execucao
	int time;
	FILE *fp;

	fp = fopen("borwein_sequencial.out", "w+");

	mpf_set_default_prec(BIT_PRECISION);

	mpf_t a0;
	mpf_t a1;
	mpf_t y0;
	mpf_t y1;
	mpf_t aux, aux2, aux3, sqrt, prec;
	unsigned long int power_of_two;
	int conv;

	mpf_init_set_ui(a0, 6);
	mpf_init(a1);
	mpf_init(y0);
	mpf_init(y1);
	mpf_init(aux);
	mpf_init(aux2);
	mpf_init(aux3);
	mpf_init(sqrt);

	printf("Calculando precisao...\n");

	// Calculo da precisao (condicao de parada)
	// prec = 1 / 10 ^ 10 000 000 
	mpf_init_set_ui(aux2, 1);
	mpf_init_set_ui(prec, 10);
	mpf_pow_ui(prec, prec, 10000000);
	mpf_ui_div(prec, 1, prec);

	gettimeofday(&start, NULL);
	printf("Calculando pi...\n");

	// Raiz(2)
	mpf_sqrt_ui(sqrt, 2);
	// aux = 4 * raiz(2)
	mpf_mul_ui(aux, sqrt, 4);
	// a0 = 6 - 4 * raiz(2)
	mpf_sub(a0, a0, aux);

	// y0 = raiz(2) - 1
	mpf_sub_ui(y0, sqrt, 1);

	int i;
	for (i = 0; conv > 0 || i == 0; i++) {
		// Calculando y1:
		// y0 = y0 ^ 4
		mpf_pow_ui(y0, y0, 4);
		// aux = 1 - y0
		mpf_ui_sub(aux, 1, y0);
		// aux = raiz(aux)
		mpf_sqrt(aux, aux);
		// aux = raiz(aux)
		mpf_sqrt(aux, aux);
		// aux2 = 1 - aux
		mpf_ui_sub(aux2, 1, aux);
		// aux3 = 1 + aux
		mpf_add_ui(aux3, aux, 1);
		// y1 = aux2 / aux3
		mpf_div(y1, aux2, aux3);

		// Calculando a1:
		// aux = y1 ^ 2
		mpf_pow_ui(aux, y1, 2);
		// aux = aux + y1
		mpf_add(aux, aux, y1);
		// aux = aux +1
		mpf_add_ui(aux, aux, 1);

		// power_of_two = 2 * i + 3
		power_of_two = pow(2, (2 * i) + 3);
		// aux2 = y1 * power_of_two
		mpf_mul_ui(aux2, y1, power_of_two);

		// aux3 = aux * aux2
		mpf_mul(aux3, aux, aux2);

		// aux = y1 + 1
		mpf_add_ui(aux, y1, 1);
		// aux = aux ^ 4
		mpf_pow_ui(aux, aux, 4);
		// aux = aux * a0
		mpf_mul(aux, aux, a0);

		// a1 = aux - aux3
		mpf_sub(a1, aux, aux3);

		// Verifica a precisao
		mpf_sub(aux, a1, a0);
		mpf_abs(aux, aux);
		conv = mpf_cmp(aux, prec);

		mpf_set(y0, y1);
		mpf_set(a0, a1);

	}

	printf("Pronto.\n");

	gettimeofday(&end, NULL);
	time = end.tv_sec - start.tv_sec;

	printf("\nNumero de iteracoes: %d\nTempo para calculo: %d segundos\n", i, time);
	printf("Escrevendo saida no arquivo...\n");
	mpf_ui_div(a0, 1, a0);
	mpf_out_str(fp, 10, 0, a0);
	printf("Pronto: borwein_sequencial.out\n");

	// Libera a memoria
	mpf_clear(a0);
	mpf_clear(a1);
	mpf_clear(y0);
	mpf_clear(y1);
	mpf_clear(aux);
	mpf_clear(aux2);
	mpf_clear(aux3);
	mpf_clear(sqrt);
	mpf_clear(prec);

	return EXIT_SUCCESS;
}
